EDA-1
Andriy Koval
May 7, 2020
Abstract
To be added
Environment
library(ggplot2) #For graphing
library(magrittr) #Pipes
library(dplyr) # for shorter function names. but still prefer dplyr:: stems
library(knitr) # dynamic documents
library(rmarkdown) # dynamic
library(kableExtra) # enhanced tables, see http://haozhu233.github.io/kableExtra/awesome_table_in_html.html
# library(TabularManifest) # exploratory data analysis, see https://github.com/Melinae/TabularManifest
requireNamespace("knitr", quietly=TRUE)
requireNamespace("scales", quietly=TRUE) #For formating values in graphs
requireNamespace("RColorBrewer", quietly=TRUE)
requireNamespace("dplyr", quietly=TRUE)
requireNamespace("DT", quietly=TRUE) # for dynamic tables
# requireNamespace("plyr", quietly=TRUE)
# requireNamespace("reshape2", quietly=TRUE) #For converting wide to long
# requireNamespace("mgcv, quietly=TRUE) #For the Generalized Additive Model that smooths the longitudinal graphs.config <- config::get()
source("./scripts/common-functions.R") # reporting functions and quick views
# source("./scripts/graphing/graph-presets.R") # font and color conventions
# source("./scripts/graphing/graph-support.R") # font and color conventionsggplot2::theme_set(ggplot2::theme_bw())compute_epi_timeline <- function(d, n_deaths_first_day = 1) { #}, d_country ){
# browser()
# d <- ds_cgrt %>%
# # filter(country_code %in% c("ITA","FRA") ) %>%
# filter(country_code %in% c("AFG") ) %>%
# select(country_code, date, n_cases, n_deaths)
#
d_out <- d %>%
# dplyr::filter(country_code %in% unique(d_country$id)) %>%
dplyr::group_by(country_code) %>%
dplyr::mutate(
# this solution might be vulnerable to cases where some intermediate dates are missed
n_deaths_cum = cumsum(tidyr::replace_na(n_deaths,0))
,n_cases_cum = cumsum(tidyr::replace_na(n_cases,0))
,cutoff_death = n_deaths_cum >= 1
,cutoff_case = n_cases_cum >= 1
,days_since_1death = cumsum(tidyr::replace_na(cutoff_death,0))
,days_since_1case = cumsum(tidyr::replace_na(cutoff_case,0))
,date_of_1death = lubridate::as_date(ifelse(days_since_1death==1,date, NA))
,date_of_1case = lubridate::as_date(ifelse(days_since_1case==1,date, NA))
,date_of_1death = min(date_of_1death, na.rm =T)
,date_of_1case = min(date_of_1case, na.rm =T)
,days_since_1death = date - date_of_1death
,days_since_1case = date - date_of_1case
,n_deaths_cum_per_1m = n_deaths_cum/n_population_2018*1000000
,n_cases_cum_per_1m = n_cases_cum/ n_population_2018*1000000
) %>%
dplyr::ungroup() %>%
# dplyr::filter(epi_timeline > 0) %>%
dplyr::mutate(
days_since_exodus = date - lubridate::date("2020-01-13") # first case outside of china
,days_since_pandemic = date - lubridate::date("2020-03-11") # WHO declares pandemic
) %>%
select(-cutoff_death, - cutoff_case, -date_of_1death, -date_of_1case)
return(d_out)
}
# for testing the function:
# d_out <- ds0 %>% filter(country_code == "ITA") %>%
# select(
# country_code, date,n_cases, n_deaths, ConfirmedDeaths, ConfirmedCases
# ) %>%
# compute_epi_timeline()Data
The data comes from [OxCGRT][] and [OECD][] databases
# list of focal countries in OECD database
ds_country <-
readr::read_csv(
config$path_country
) %>%
dplyr::filter(desired)
#
# # ECDC
# # path_save <- paste0("./data-unshared/derived/ocdc-",Sys.Date(),".csv")
ds_covid <- readr::read_csv(config$path_input_covid,)
# ds_covid %>% glimpse()
ds_country_codes <- readr::read_csv(config$path_country_codes)
# OECD
file_path <- list.files(config$path_oecd_clean,full.names = T,recursive = T,pattern = ".rds$")
dto <- list()
for(i in seq_along(file_path)){
file_name <- basename(file_path[i]) %>% stringr::str_replace(".rds","")
dto[[file_name]] <- readr::read_rds(file_path[i])
}
# str(dto,max.level = 1)
ls_health_resources <- dto$health_resources
# ls_health_resources %>% str(1)
ds_hr <- ls_health_resources$data_agg
# ds_hr %>% glimpse()
# OxCGRT
ds_cgrt <- readr::read_rds("./data-unshared/derived/OxCGRT.rds")
# ds_cgrt %>% glimpse()
# n_distinct(ds_cgrt$country_code)Data Tweaks
ds0 <- ds_covid %>%
compute_epi_timeline() %>%
dplyr::left_join(
ds_cgrt
,by = c("date", "country_code")
) %>%
dplyr::left_join(
ds_country_codes,
by = c("country_code" = "country_code3")
)
d_out <- ds0 %>% filter(country_code == "ITA")Q1 - Shape of the trend
How does COVID-19 progress in each country?
# ds0 %>% glimpse()
d1 <- ds0 %>%
filter(country_code %in% ds_country$id)
g1 <- d1 %>%
ggplot(aes(
x = days_since_exodus
,y = n_cases_cum
# ,y =n_cases_cum_per_1m
# ,y = n_deaths_cum
# ,y = n_deaths_cum_per_1m
))+
geom_line()+
# geom_line(aes(y=StringencyIndex), color = "red")+
facet_wrap(~country_label, scale = "free")+
geom_point(data = d1 %>% filter(days_since_1case == 1), size = 2, fill = "#1b9e77", color = "black", alpha = .5, shape = 21)+
geom_point(data = d1 %>% filter(days_since_1death == 1), size = 2, fill = "#d95f02", color = "black", alpha = .5, shape = 21)+
labs(
title = "Timeline of COVID-19 "
, y = "Cumulative Cases", x = "Days since first case outside of China (Jan 13, 2020)"
, caption = "(first dot) = 1st confirmed case, (second dot) = 1st confirmed death, (dashed line) = pandemic announced by WHO"
)+
geom_vline(xintercept = 58, linetype = "dashed")
cat("\n## Cases\n")Cases
g1
cat("\n## Cases per 1m\n")Cases per 1m
g1 + aes(y = n_cases_cum_per_1m)+labs(y = "Cumulative Cases per 1 mil")
cat("\n## Deaths\n")Deaths
g1 + aes(y = n_deaths_cum)+labs(y = "Cumulative Deaths")
cat("\n## Deaths per 1m\n")Deaths per 1m
g1 + aes(y = n_deaths_cum_per_1m)+labs(y = "Cumulative Deaths per 1 mil")
Q1 - Within Country
focal_vars <- c( "n_cases_cum", "n_cases_cum_per_1m", "n_deaths_cum", "n_deaths_cum_per_1m",
"StringencyIndex")
ds1 <- ds0 %>%
filter(country_code %in% ds_country$id) %>%
# dplyr::filter(country_code %in% c("ITA","FRA")) %>%
dplyr::select(country_code, country_label, days_since_exodus, days_since_1case,
days_since_1death,
n_cases_cum, n_cases_cum_per_1m, n_deaths_cum, n_deaths_cum_per_1m,
StringencyIndex
) %>%
tidyr::pivot_longer(cols = focal_vars, values_to = "value", names_to = "metric")
print_one_wrap <- function(d, country = "ITA"){
# d <- ds1; country = "ITA"
d1 <- d %>% filter(country_code == country)
g1 <- d1 %>%
ggplot(aes(x = days_since_exodus, y = value))+
geom_line()+
geom_point(
data = d1 %>% filter(days_since_1case == 1),
size = 2, fill = "#1b9e77", color = "black", alpha = .5, shape = 21)+
geom_point(
data = d1 %>% filter(days_since_1death == 1),
size = 2, fill = "#d95f02", color = "black", alpha = .5, shape = 21)+
geom_vline(xintercept = 58, linetype = "dashed")+
facet_wrap(country_label ~ metric, scale = "free_y",ncol = 5)
g1
}
# ds1 %>% print_one_wrap(country = 'ITA')
countries <- unique(ds1$country_code)
for(country_i in countries){
cat("\n## ", country_i,"\n")
ds1 %>% print_one_wrap(country = country_i) %>% print()
cat("\n")
}AUS
AUT
BEL
CAN
CHE
CZE
DEU
DNK
ESP
EST
FIN
FRA
GBR
GRC
HUN
IRL
ISL
ISR
ITA
LTU
LVA
NLD
NOR
NZL
POL
PRT
RUS
SVK
SVN
SWE
TUR
USA
ZAF
Q2 - Response Stringency
What was the trend of the response to COVID-10 by each country?
# What the trend response to COVID-10 by each country?
d1 <- ds0 %>%
filter(country_code %in% ds_country$id)
g1 <- d1 %>%
ggplot(aes(x = days_since_exodus, y = StringencyIndex))+
geom_line()+
geom_point(data = d1 %>% filter(days_since_1case == 1), size = 2, fill = "#1b9e77", color = "black", alpha = .5, shape = 21)+
geom_point(data = d1 %>% filter(days_since_1death == 1), size = 2, fill = "#d95f02", color = "black", alpha = .5, shape = 21)+
facet_wrap(~country_label)+
labs(
title = "Timeline of OECD countries' respones to COVID-19 as measured by the Stringency Index"
,y = "Stringency Index", x = "Days since first case outside of China (Jan 13, 2020)"
, caption = "First dot = 1st confired case, Second dot = 1st confirmed death, line = Pandemic announced by WHO"
)+
geom_vline(xintercept = 58, linetype = "dashed")
g1
d2 <- ds0 %>%
filter(country_code %in% ds_country$id)
g2 <- d2 %>%
filter(country_code %in% ds_country$id) %>%
# filter(country_code == "ITA") %>%
ggplot(aes(x = days_since_exodus, y = StringencyIndex, group = country_label))+
geom_line( alpha = .2)+
geom_point(data = d2 %>% filter(days_since_1case == 1), size = 2, fill = "#1b9e77", color = "black", alpha = .5, shape = 21)+
geom_point(data = d2 %>% filter(days_since_1death == 1), size = 2, fill = "#d95f02", color = "black", alpha = .5, shape = 21)+
labs(
title = "Timeline of OECD countries' respones to COVID-19 as measured by the Stringency Index"
,y = "Stringency Index", x = "Days since first case outside of China (Jan 13, 2020)"
)+
geom_vline(xintercept = 58, linetype = "dashed")
g2 <- plotly::ggplotly(g2)
g2Q3 - Response metrics
# d1 <- ds0 %>%
d_out <- ds0 %>%
filter(country_code == "ITA") %>%
select(country_code, date,n_cases_cum, n_deaths_cum, days_since_1case, days_since_1death)
# Deaths 30 days after 1st death
d1 <- ds0 %>%
group_by(country_code) %>%
dplyr::filter(days_since_1death == 30) %>%
dplyr::select(country_code, n_deaths_cum, n_population_2018) %>%
dplyr::rename(n_deaths_30days_since_1death = n_deaths_cum) %>%
dplyr::mutate(n_deaths_30days_since_1death_per100k = n_deaths_30days_since_1death/n_population_2018*100000 ) %>%
dplyr::select(-n_population_2018)
# Deaths 60 days after 1st death
d2 <- ds0 %>%
group_by(country_code) %>%
dplyr::filter(days_since_1death == 60) %>%
dplyr::select(country_code, n_deaths_cum,n_population_2018) %>%
dplyr::rename(n_deaths_60days_since_1death = n_deaths_cum) %>%
dplyr::mutate(
n_deaths_60days_since_1death_per100k = n_deaths_60days_since_1death/n_population_2018*100000
) %>%
dplyr::select(-n_population_2018)
# Cases 30 days after 1st case
d3 <- ds0 %>%
group_by(country_code) %>%
dplyr::filter(days_since_1case == 30) %>%
dplyr::select(country_code, n_cases_cum, n_population_2018) %>%
dplyr::rename(n_cases_30days_since_1case = n_cases_cum) %>%
dplyr::mutate(
n_cases_30days_since_1case_per100k = n_cases_30days_since_1case/n_population_2018*100000
) %>%
dplyr::select(-n_population_2018)
# Cases 60 days after 1st case
d4 <- ds0 %>%
group_by(country_code) %>%
dplyr::filter(days_since_1case == 60) %>%
dplyr::select(country_code, n_cases_cum, n_population_2018) %>%
dplyr::rename(n_cases_60days_since_1case = n_cases_cum) %>%
dplyr::mutate(
n_cases_60days_since_1case_per100k = n_cases_60days_since_1case/n_population_2018*100000
) %>%
dplyr::select(-n_population_2018)
# ds0 %>% filter(country_ == "LVA")
ds_response <- list(d1,d2,d3,d4) %>% Reduce(function(a,b) dplyr::full_join(a,b), .)
ds_response <- ds_response %>%
dplyr::left_join(
ds_covid %>% distinct(country_code, geo_id, country)
) %>%
dplyr::filter(!is.na(country_code))
ds_response %>% neat_DTsession information
For the sake of documentation and reproducibility, the current report was rendered in the following environment. Click the line below to expand.
Environment
- Session info -------------------------------------------------------------------------------------------------------
setting value
version R version 3.6.3 (2020-02-29)
os Windows 10 x64
system x86_64, mingw32
ui RTerm
language (EN)
collate English_United States.1252
ctype English_United States.1252
tz America/New_York
date 2020-05-08
- Packages -----------------------------------------------------------------------------------------------------------
package * version date lib source
assertthat 0.2.1 2019-03-21 [1] CRAN (R 3.6.2)
backports 1.1.5 2019-10-02 [1] CRAN (R 3.6.1)
callr 3.4.3 2020-03-28 [1] CRAN (R 3.6.3)
cli 2.0.2 2020-02-28 [1] CRAN (R 3.6.3)
codetools 0.2-16 2018-12-24 [2] CRAN (R 3.6.3)
colorspace 1.4-1 2019-03-18 [1] CRAN (R 3.6.1)
config 0.3 2018-03-27 [1] CRAN (R 3.6.3)
crayon 1.3.4 2017-09-16 [1] CRAN (R 3.6.2)
crosstalk 1.0.0 2016-12-21 [1] CRAN (R 3.6.2)
data.table 1.12.8 2019-12-09 [1] CRAN (R 3.6.2)
desc 1.2.0 2018-05-01 [1] CRAN (R 3.6.2)
devtools 2.3.0 2020-04-10 [1] CRAN (R 3.6.3)
digest 0.6.25 2020-02-23 [1] CRAN (R 3.6.3)
dplyr * 0.8.5 2020-03-07 [1] CRAN (R 3.6.3)
DT 0.13 2020-03-23 [1] CRAN (R 3.6.3)
ellipsis 0.3.0 2019-09-20 [1] CRAN (R 3.6.2)
evaluate 0.14 2019-05-28 [1] CRAN (R 3.6.2)
fansi 0.4.1 2020-01-08 [1] CRAN (R 3.6.2)
farver 2.0.3 2020-01-16 [1] CRAN (R 3.6.2)
fastmap 1.0.1 2019-10-08 [1] CRAN (R 3.6.2)
fs 1.3.1 2019-05-06 [1] CRAN (R 3.6.2)
generics 0.0.2 2018-11-29 [1] CRAN (R 3.6.2)
ggplot2 * 3.2.1 2019-08-10 [1] CRAN (R 3.6.2)
glue 1.4.0 2020-04-03 [1] CRAN (R 3.6.3)
gtable 0.3.0 2019-03-25 [1] CRAN (R 3.6.2)
hms 0.5.3 2020-01-08 [1] CRAN (R 3.6.2)
htmltools 0.4.0 2019-10-04 [1] CRAN (R 3.6.2)
htmlwidgets 1.5.1 2019-10-08 [1] CRAN (R 3.6.2)
httpuv 1.5.2 2019-09-11 [1] CRAN (R 3.6.2)
httr 1.4.1 2019-08-05 [1] CRAN (R 3.6.2)
jsonlite 1.6.1 2020-02-02 [1] CRAN (R 3.6.2)
kableExtra * 1.1.0 2019-03-16 [1] CRAN (R 3.6.3)
knitr * 1.28 2020-02-06 [1] CRAN (R 3.6.2)
labeling 0.3 2014-08-23 [1] CRAN (R 3.6.0)
later 1.0.0 2019-10-04 [1] CRAN (R 3.6.2)
lazyeval 0.2.2 2019-03-15 [1] CRAN (R 3.6.2)
lifecycle 0.2.0 2020-03-06 [1] CRAN (R 3.6.3)
lubridate 1.7.8 2020-04-06 [1] CRAN (R 3.6.3)
magrittr * 1.5 2014-11-22 [1] CRAN (R 3.6.2)
memoise 1.1.0 2017-04-21 [1] CRAN (R 3.6.2)
mime 0.9 2020-02-04 [1] CRAN (R 3.6.2)
munsell 0.5.0 2018-06-12 [1] CRAN (R 3.6.2)
pillar 1.4.3 2019-12-20 [1] CRAN (R 3.6.2)
pkgbuild 1.0.6 2019-10-09 [1] CRAN (R 3.6.2)
pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 3.6.2)
pkgload 1.0.2 2018-10-29 [1] CRAN (R 3.6.2)
plotly 4.9.2 2020-02-12 [1] CRAN (R 3.6.2)
prettyunits 1.1.1 2020-01-24 [1] CRAN (R 3.6.2)
processx 3.4.2 2020-02-09 [1] CRAN (R 3.6.2)
promises 1.1.0 2019-10-04 [1] CRAN (R 3.6.2)
ps 1.3.2 2020-02-13 [1] CRAN (R 3.6.2)
purrr 0.3.4 2020-04-17 [1] CRAN (R 3.6.3)
R6 2.4.1 2019-11-12 [1] CRAN (R 3.6.2)
RColorBrewer 1.1-2 2014-12-07 [1] CRAN (R 3.6.0)
Rcpp 1.0.4.6 2020-04-09 [1] CRAN (R 3.6.3)
readr 1.3.1 2018-12-21 [1] CRAN (R 3.6.2)
remotes 2.1.1 2020-02-15 [1] CRAN (R 3.6.2)
rlang 0.4.5 2020-03-01 [1] CRAN (R 3.6.3)
rmarkdown * 2.1 2020-01-20 [1] CRAN (R 3.6.2)
rprojroot 1.3-2 2018-01-03 [1] CRAN (R 3.6.2)
rstudioapi 0.11 2020-02-07 [1] CRAN (R 3.6.2)
rvest 0.3.5 2019-11-08 [1] CRAN (R 3.6.2)
scales 1.1.0 2019-11-18 [1] CRAN (R 3.6.2)
sessioninfo 1.1.1 2018-11-05 [1] CRAN (R 3.6.2)
shiny 1.4.0 2019-10-10 [1] CRAN (R 3.6.2)
stringi 1.4.6 2020-02-17 [1] CRAN (R 3.6.2)
stringr 1.4.0 2019-02-10 [1] CRAN (R 3.6.2)
testthat 2.3.2 2020-03-02 [1] CRAN (R 3.6.3)
tibble 3.0.1 2020-04-20 [1] CRAN (R 3.6.3)
tidyr 1.0.2 2020-01-24 [1] CRAN (R 3.6.2)
tidyselect 1.0.0 2020-01-27 [1] CRAN (R 3.6.2)
usethis 1.6.0 2020-04-09 [1] CRAN (R 3.6.3)
vctrs 0.2.4 2020-03-10 [1] CRAN (R 3.6.3)
viridisLite 0.3.0 2018-02-01 [1] CRAN (R 3.6.2)
webshot 0.5.2 2019-11-22 [1] CRAN (R 3.6.3)
withr 2.1.2 2018-03-15 [1] CRAN (R 3.6.2)
xfun 0.12 2020-01-13 [1] CRAN (R 3.6.2)
xml2 1.2.2 2019-08-09 [1] CRAN (R 3.6.2)
xtable 1.8-4 2019-04-21 [1] CRAN (R 3.6.2)
yaml 2.2.1 2020-02-01 [1] CRAN (R 3.6.2)
[1] C:/Users/an499583/Documents/R/win-library/3.6
[2] C:/Users/an499583/Documents/R/R-3.6.3/library